/*
 * To change this license header, choose License Headers in Project Properties.
 * To change this template file, choose Tools | Templates
 * and open the template in the editor.
 */
package com.bioinfseq;

import java.io.BufferedReader;
import java.io.FileNotFoundException;
import java.io.FileReader;
import java.io.IOException;
import java.util.logging.Level;
import java.util.logging.Logger;
import org.biojava3.alignment.Alignments;
import org.biojava3.alignment.Alignments.PairwiseSequenceAlignerType;
import org.biojava3.alignment.SimpleGapPenalty;
import org.biojava3.alignment.SimpleSubstitutionMatrix;
import org.biojava3.alignment.template.PairwiseSequenceAligner;
import org.biojava3.alignment.template.SequencePair;
import org.biojava3.core.sequence.DNASequence;
import org.biojava3.core.sequence.RNASequence;
import org.biojava3.core.sequence.compound.AmbiguityDNACompoundSet;
import org.biojava3.core.sequence.compound.AminoAcidCompound;
import org.biojava3.core.sequence.compound.AminoAcidCompoundSet;
import org.biojava3.core.sequence.compound.DNACompoundSet;
import org.biojava3.core.sequence.compound.NucleotideCompound;
import org.biojava3.core.sequence.template.Sequence;
import org.biojava3.core.sequence.transcription.TranscriptionEngine;

/**
 *
 * @author Bartek
 */
public class Application {

    public static void main(String[] args) {
        
        String inputPath;
        if (args.length > 0) {
            inputPath = args[0];
        }
        else {
            inputPath = "dna_input.txt";
        }
        
        BufferedReader reader = null;
        String firstSequence = "", secondSequence = "", option = "", penalties = "";
        short openPenalty = 0, extensionPenalty = 0;
        SimpleSubstitutionMatrix<NucleotideCompound> substitutionMatrix = null;
        SimpleSubstitutionMatrix<AminoAcidCompound> substitutionMatrixAmino = null;
        try {
            reader = new BufferedReader(new FileReader(inputPath));
            option = reader.readLine();
            penalties = reader.readLine();
            firstSequence = reader.readLine();
            secondSequence = reader.readLine();
            if ("DNA".equals(option)) {
                substitutionMatrix
                        = new SimpleSubstitutionMatrix<>(DNACompoundSet.getDNACompoundSet(), reader, "mm");
            }
            if ("RNA".equals(option)) {
                substitutionMatrixAmino
                        = new SimpleSubstitutionMatrix<>(AminoAcidCompoundSet.getAminoAcidCompoundSet(), reader, "mm");
            }

        } catch (FileNotFoundException ex) {
            Logger.getLogger(Application.class.getName()).log(Level.SEVERE, null, ex);
        } catch (IOException ex) {
            Logger.getLogger(Application.class.getName()).log(Level.SEVERE, null, ex);
        } finally {
            try {
                if (reader != null) {
                    reader.close();
                }
            } catch (IOException ex) {
                Logger.getLogger(Application.class.getName()).log(Level.SEVERE, null, ex);
            }
        }
        
        if ("DNA".equals(option)) {
            DNASequence target = new DNASequence(firstSequence,
                    AmbiguityDNACompoundSet.getDNACompoundSet());

            DNASequence query = new DNASequence(secondSequence,
                    AmbiguityDNACompoundSet.getDNACompoundSet());

            CustomPenalty cPenalty = new CustomPenalty(openPenalty, extensionPenalty);

            PairwiseSequenceAligner<DNASequence, NucleotideCompound> localAligner
                    = Alignments.getPairwiseAligner(query, target,
                            PairwiseSequenceAlignerType.LOCAL, cPenalty, substitutionMatrix);

            PairwiseSequenceAligner<DNASequence, NucleotideCompound> globalAligner
                    = Alignments.getPairwiseAligner(query, target,
                            PairwiseSequenceAlignerType.GLOBAL, cPenalty, substitutionMatrix);

            printAlignment("Local alignment", localAligner);
            printAlignment("Global alignment", globalAligner);

        } else if ("RNA".equals(option)) {
            RNASequence target = new RNASequence(firstSequence);
            RNASequence query = new RNASequence(secondSequence);
            CustomPenalty cPenalty = new CustomPenalty(openPenalty, extensionPenalty);

            TranscriptionEngine engine = TranscriptionEngine.getDefault();
            Sequence<AminoAcidCompound> firstProtein = engine.getRnaAminoAcidTranslator().createSequence(target);
            Sequence<AminoAcidCompound> secondProtein = engine.getRnaAminoAcidTranslator().createSequence(query);

            System.out.println("First protein:");
            System.out.println(firstProtein.toString());
            System.out.println("Second protein:");
            System.out.println(secondProtein.toString());
            //SubstitutionMatrix<AminoAcidCompound> aminoMatrix = SubstitutionMatrixHelper.getBlosum65();

            PairwiseSequenceAligner<Sequence<AminoAcidCompound>, AminoAcidCompound> localAlignment
                    = Alignments.getPairwiseAligner(firstProtein, secondProtein, PairwiseSequenceAlignerType.LOCAL, cPenalty, substitutionMatrixAmino);
            
            SequencePair<Sequence<AminoAcidCompound>, AminoAcidCompound> pairLocal = localAlignment.getPair();

            PairwiseSequenceAligner<Sequence<AminoAcidCompound>, AminoAcidCompound> globalAlignment
                    = Alignments.getPairwiseAligner(firstProtein, secondProtein, PairwiseSequenceAlignerType.GLOBAL, cPenalty, substitutionMatrixAmino);

            SequencePair<Sequence<AminoAcidCompound>, AminoAcidCompound> pairGlobal = globalAlignment.getPair();

            System.out.println("Local alignment:");
            System.out.println(pairLocal.toString());
            System.out.println("Global aligment:");
            System.out.println(pairGlobal.toString());

            firstProtein = pairLocal.getQuery();
            secondProtein = pairLocal.getTarget();

            System.out.println("First protein:");
            System.out.println(firstProtein.toString());
            System.out.println("Second protein:");
            System.out.println(secondProtein.toString());

        }

    }
    
    private static void printAlignment(String name, PairwiseSequenceAligner<DNASequence, NucleotideCompound> aligner) {
        System.out.println(name);
        System.out.println(aligner.getPair());
        System.out.println("Edit distance: " + aligner.getDistance());
        System.out.println("Similarity: " + aligner.getSimilarity());
    }

}
